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IMPROVED DETERMINATION OF VECTOR LITHOSPHERIC 
MAGNETIC ANOMALIES FROM MAGSAT DATA 


This final report describes scientific contributions made under NASA 
grant NAG5-1972 in developing new methods to isolate and map vector 
magnetic anomalies from measurements made by Magsat. In addition to 
the objective of the proposal, the isolation and mapping of equatorial 
vector lithospheric Magsat anomalies, we have tackled isolation of polar 
ionospheric fields during the award period. Significant progress has also 
been made in isolation of polar AZ component and scalar anomalies as well 
as integration and synthesis of various techniques of removing equatorial 
and polar ionospheric effects. This research has resulted in one invited 
published abstract and one manuscript is under preparation discussing the 
outcome of the research (to be submitted to a reviewed journal). 

The significant contributions of this research are: 1) development of 
empirical/analytical techniques in modeling ionospheric fields in Magsat 
data and their removal from uncorrected anomalies to obtain better 
estimates of lithospheric anomalies. This task has been accomplished for 
equatorial AX, AZ, and AB component and polar AZ and AB component 
measurements; 2) integration of important processing techniques 
developed during the last decade with the newly developed technologies of 
ionospheric field modeling into an optimum processing scheme, and 3) 
implementation of the above processing scheme to map the most robust 
magnetic anomalies of the lithosphere (components as well as scalar). 

The results of the research have been documented in an invited 
paper presented at the Fall 1992 American Geophysical Union meeting 
under: 

Ravat, D., R. Langel, M. Purucker, T. Sabaka, J. Arkani-Hamed, D. 
Alsdorf, A new approach for isolation of lithospheric magnetic 
anomalies from Magsat, Invited paper, Eos Trans. AGU, 73, 
p.140, 1992. 



These results are also being documented, to be published in a reviewed 
journal (under preparation), as: 

Ravat, D., R. Langel, M. Purucker, J. Arkani-Hamed, and D. Alsdorf, 
Vector and scalar Magsat magnetic anomaly maps for geologic 
interpretation, in preparation. 

The text of this manuscript describes the details of the procedures 
developed during the study and their results and, hence, is included as the 
remainder of the final report. 


VECTOR AND SCALAR MAGS AT MAGNETIC ANOMALY MAPS FOR GEOLOGIC 

INTERPRETATION 


INTRODUCTION 

Satellite-derived magnetic anomalies are proving useful in mapping 
and understanding geologic evolution of large-scale tectonic features on 
the earth. These data have been used to map large-scale continental 
tectonic provinces (Frey, 1982; Arkani-Hamed and Strangway, 1985a; 
Ravat, 1989; Toft et al., 1992) and mineral resources (Taylor et al., 1992; 
Ravat et al., 1993), oceanic uplifts (Fullerton et al., 1989; Antoine and 
Moyes, 1992), magnetic quiet zones on the seafloor (LaBrecque and 
Raymond, 1985), and the nature and evolution of the Mesozoic break-up of 
Pangea (Frey et al., 1983; Galdeano, 1983; von Frese et al., 1986; Ravat et 
al., 1992). 

Most researchers working with satellite magnetic data have used 
total intensity magnetic anomalies from the POGO and Magsat missions 
because of the relatively low error envelope of these data when compared 
with vector component data. Although Magsat measured three-orthogonal- 
component vector magnetic fields at an average altitude of approximately 
400 km, the accuracy of the vector measurements depends largely on how 
well the directions of the magnetometer sensors are known (2-20 arc-sec). 
Thus, vector observations contain an inherent source of error (0.5-5 nT) 
beyond those in total intensity data (1-2 nT). The problem of accurate 
attitude determination is exacerbated by discontinuities or jumps in 
attitude solutions that result from changing the configuration of the 
attitude sensing instruments many times an orbit (Mayhew et al., 1985). 
Moreover, geologic anomalies are the residual of various large and small 
magnetic contributions to the measurements (e.g.. Earth's main magnetic 
field, effects of magnetospheric and ionospheric currents, etc.). The 
precision of magnetic fields which represent geologic sources (referred 
interchangeably to as 'geologic', 'crustal', or 'lithospheric' magnetic 
anomalies in this paper) is thus limited by how well these other magnetic 
fields can be removed. 

The problem faced is illustrated by the maps of Figures 1 and 2. As 


described by Langel et al. (1982), the Magsat satellite acquired data only 
at twilight local times. Thus a natural division of that data into 
independent . subsets is to treat dawn and dusk data separately (Arkani- 
Hamed et al., 1985; Yanagisawa and Kono, 1985). Figures 1 and 2 are 
residual maps of the total intensity anomalies of Magsat data at dawn and 
at dusk. The figures show ionospheric contamination in the form of 
anomaly bands that follow magnetic and geomagnetic latitudes in the low 
and high latitude regions, respectively. Data utilized in these maps has 
been selected from magnetically quiet times and low pass filtered 
(procedures to be described in detail in a later paragraph). A similar 
representation of the three vector component maps shows interesting 
patterns (not shown here). The AY map shows distinct positive residuals to 
the north and negative residuals to the south of the dip equator; similarly, 
the AZ map shows negative residuals to the north and positive residuals to 
the south of the dip equator; and the AX map shows a negative residual 
along the dip equator. Maeda et al. (1982, 1985) showed that this 
distinctive AY pattern is due to a meridional current associated with the 
equatorial electrojet. The AX and AZ patterns are consistent with an 
eastward flowing equatorial electrojet below the satellite. Data from dawn, 
especially AY, are relatively free from these effects, except as introduced 
by errors in deriving the residuals by removing a main field model (Langel 
et al., 1993). 

Fortunately, many of the problems mentioned and illustrated in 
Figures 1 and 2 are not insurmountable. During the past decade, an arsenal 
of techniques has been developed and applied to both the scalar and 
vector data sets in attempts to isolate the geologic anomalies (Yanagisawa 
and Kono, 1984, 1985; Arkani-Hamed and Strangway, 1985a,b, 1986; 
Nakagawa and Yukutake, 1984, 1985; Nakagawa et al., 1985; Taylor and 
Frawley, 1987; von Frese et al., 1988; Cohen, 1989; Ravat, 1989; Cohen and 
Achache, 1990; Alsdorf, 1991; Purucker, 1991; Alsdorf et al., 1993; Langel 
et al., 1993; Ravat and Hinze, 1993). The authors of the present paper have 
separately been participants in developing portions of that arsenal of 
techniques. Recently, Langel (1992) called for a procedure in which several 
of these techniques are applied in a logical order. Partially in response to 
that suggestion, during the summer of 1992 the authors gathered at 
Goddard Space Flight Center in order to work as a team. The intent of the 
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present paper is to present the results of a concentrated team effort to 
collect together, refine, and, where necessary, extend existing techniques 
with the goal of deriving the best possible scalar and, especially, vector 
anomaly maps from Magsat data under the constraint of having to 
complete most of the work in a short period of time. The maps presented 
are not likely the ultimate vector anomaly maps: further improvements in 
isolation and resolution may be possible. Indeed, during this project we 
have identified areas of research and techniques (mainly for improved 
modeling of ionospheric and magnetospheric fields) that could possibly 
yield better resolved anomaly maps in the future. However, based on our 
combined experience, the maps presented here are a considerable 
improvement over those in the past and most probably describe fields 
geologic in origin which should be amenable to geologic interpretation. 

In the following sections the methodology of how these techniques 
were applied to the data will be outlined, results from various stages of the 
process will be illustrated, and the final maps presented. 


GENERAL DATA REDUCTION STRATEGIES 


General Comments 

If B(r,t) is the measured magnetic field at location r and time t then 
B(r,t) = M(r,t) + A(r) + D(r,t) + e ( 1 ) 

where M(r,t) is the field from Earth's core, the main field, A(r) is the field 
from Earth's crust, D(r,t) is the field from magnetospheric and ionospheric 
sources, including portions induced in Earth, and e is the measurement 
error. In practice, analysis of both A(r) and D(r,t) begins by removing an 
estimate of M(r,t) from B(r,t). That estimate is from some model of the 
main field. The difference, or residual field is then 


\ 


AB(r,t) = B(r,t) - M'(r,t) = A(r) + D(r,t) + q 


( 2 ) 


where M'(r,t) is an estimate of M (r,t) and q is a combination of e and 
model error. For crustal anomaly studies, D(r,t) and q constitute "noise". 

The scalar residual field, or AB, is by definition 

AB(r,t) = |B(r,t)l - lM(r,t)l. (3) 

Now, ignoring e, 

|AB| = [(AB + M)*(AB+M)]i/2 = [AB«AB+2AB«M+M*M] 1/2 . (4) 

But AB*AB in (4) is much smaller than the other terms and can be 
neglected. Furthermore, the square root can be well approximated by the 
first two terms of its Taylor expansion. As a result 

AB = AB*M/|M| (5) 

is a good approximation (Note: AB * IABI). 

Figure 3 is a flow-chart to summarize the data processing sequence. 
Important parameters used in application of a particular technique are 
given in Table 1. The following sections detail the process in the flow- 
chart and the reasons behind the choice of parameters. Choosing one 
technique over another designed for the same purpose is somewhat 
arbitrary and does not necessarily imply superiority of one technique over 
the other. Our choice is largely dictated by a judgement regarding the most 
suitable technique for the task in hand, but also involves complex 
considerations such as past experience with the technique, personal biases, 
and the democratic process. 

During this project, known results of ionospheric and magnetospheric 
physics have been utilized as much as possible. In some instances, 
physical/empirical corrections have been directly applied (e.g., corrections 
for Equatorial Electrojet fields), in others they have been applied in an 
average sense (e.g., corrections applied for fields of ionospheric origin in 
polar regions). In other instances, the cause behind a particular variation 


simply could not be unambiguously deduced (e.g., some low latitude 
seasonal-temporal effects at dawn local-time). However, even when a 
good measure of success was achieved in correcting for D(r,t), this study 
also has had to rely on statistical/signal processing methods to reduce 
discrepancies between maps. 


Selection of Quiet Data 


When deriving maps of crustal anomalies, A(r), it is sensible to begin 
with data from times of relatively magnetically quiet conditions, called 
"quiet data". How quiet depends on availability of sufficient data values 
with the chosen quiet criterion to make a map. Appropriate magnetic 
indices (see, e.g., Mayaud, 1980; Rangarajan, 1989) to choose quiet data are 
Kp (for equatorial and mid-latitude data) and AE (for polar data). 
Alternately or concurrently, one can use the variance of AB along a pass to 

judge the magnetic quietness of that pass. The premise in using the 
variance criterion is that A(r) is stationary and always present whereas 
external fields, D(r,t), are dynamic, with significant variation in amplitude 
from pass-to-pass, and always add energy (area under a residual profile 
curve) to magnetic observations. Profiles in the same location with less 
anomaly variance are least likely to be contaminated by D(r,t). The criteria 
used for the selection of quiet data are listed in Table 1. In choosing the 
variance cut-off for polar data, passes over a particular region are 
screened to ascertain that the cause of high variance was indeed D(r,t) and 
not A (r). One can determine that a pass has a large amount of D(r,t) 
contribution if, over the same track location, one finds another pass with 
significantly lower variance. 


Main Field Removal 

Accuracy in modeling M (r,t), including its secular variation is 
becoming increasingly important to the improvement of lithospheric 
component maps and to the represention of ionospheric contributions to 


magnetic observations. Improper main field representation can introduce 
systematic variations in data, i.e. the tj, that can be aliased as ionospheric 
current effects or lithospheric anomalies (Langel et al., 1993; Ravat and 
Hinze, 1993). Inclusion of even small amounts of the quietest equatorial 
dusk data in GSFC(12/83) field model (Langel and Estes, 1985a) led to a 
very small contamination (on the order of 2-5 nT) of that model in low 
latitudes in such a way as to introduce a systematic deviation in dawn 
vector Magsat data (Langel et al., 1993). For this reason, for the 
equatorial/mid-latitude segment of the data, a main field model was used 
which was derived by excluding dusk Magsat data (Dawn(6/83-6), Langel 
and Estes, 1985b). For polar data, however, the GSFC(12/83) field model 
(Langel and Estes, 1985a) was used. 

Data were despiked for elimination of spurious values; over the 
poles, all selected passes were visually examined for spurious activity. 


Removal of Vector Discontinuities 

After the removal of M'(r,t), residual data were corrected for 
attitude solution jumps (Mayhew et al., 1985) which arise in processing of 
Magsat data. The attitude determination system on Magsat incorporated 
three sensors: two star trackers and a precision sun sensor. Finding the 
spacecraft attitude at any particular time involved a process of fitting the 
data from whichever sensors were contributing data at that time. If data 
were available from any two of the three sensors, an attitude solution was 
possible. However, the alignment between the three sensors was not 
perfectly known. As a result, whenever the available combination of 
sensors changed, a small discontinuity was introduced into the attitude 
solution. Many of the resulting field discontinuities are small, within 1 - 3 
nT. A significant number, however, are larger and affect attempts to derive 
anomaly maps. We have used an empirical procedure developed by 
Purucker (1991) to adjust discontinuities in the magnetic component data. 
Illustration of the adjustments is shown in Figure 4. 


Long-wavelength Filtering 

Long-wavelength filtering is applied in this study to reduce the 
long-wavelength effects of ionospheric and magnetospheric currents. Ring- 
current is an equatorial sheet current that persists from two-to-three 
Earth radii out to Earth's magnetopause or magnetotail. Langel and 
Sweeney (1971) derived a potential function representation for the 
magnetic effect of this current that removes its major contribution in 
satellite magnetic observations. However, experience in applying this 
correction to satellite data has shown that the correction does not always 
center on the broad-scale minimum indicated by observations. As a result, 
the correction leaves behind residuals of both long- and short- 
wavelengths. Also, it is now apparent that many of the long-wavelength 
inconsistencies of equatorial Magsat data are related to equatorial 
electrojet and/or Sq fields, rather than magnetospheric fields. High-pass 
filtering, with wavelength-cut from 3000 to 5500 km according to the 
correlative nature of anomaly signal, has been implemented (Ridgway and 
Hinze, 1986; Baldwin and Frey, 1991) to remove long-wavelength residual 
artifacts from the ring-current and quiet day ionospheric currents, the Sq. 
As will be discussed, Langel et al. (1993) have successfully modeled many 
features of the magnetic field from the equatorial electrojet current 
system. This modeling required extension of the filter cut-off from about 
4000 km to 12000 km for the equatorial/mid-latitude data. The filter 
selected was a "Kaiser" filter (Kaiser, 1974). The 12000 km wavelength 
cutoff does effectively remove magnetic effects of the ring-current and Sq. 

In polar regions, however, modeling of D(r,t) is much more 
complicated and still in its infancy (for the purpose of deriving lithospheric 
magnetic anomalies). Hence, these data are filtered with a 4000 km 
wavelength high- pass Kaiser filter (Figure 3, Table 1). Similarly, scalar 
equatorial data have been filtered with the 4000 km cutoff because a 
different approach has been used to remove ionospheric effects in these 
data. 


Identification of common signals: Correlation filtering 


At several points in the procedure adopted for isolating A(r) a 
technique for identifying the common features in two data sets is required. 
One instance is to find the common signal from two satellite passes along 
identical, or nearly identical, tracks. The other is to find the common 
signal in two independent residual maps, either of some component or 
AB(r,t) or of A(r). 

Pass-bv-pass correlation was used successfully by Alsdorf (1991) who 
pointed out that given the Fourier analyses of two signals of equal length, 
e.g., 

1 N - ! 

S 1 = N X A k [sin(kx) + cos(kx) ] , 

k=0 

! N-l 

s 2 = n Z B k [sin(kx + t k ) + cos(kx + tk) ] , (6) 

k=0 

then the correlation coefficient between the k 1 * 1 harmonics of the two 
analyses is given by 

pk = cos(tk) (7) 

independently of the amplitudes of the Fourier components. 

Spherical harmonic covariant correlation was developed by Arkani-Hamed 
and Strangway (1986) for analysis of scalar anomaly data but is applicable 
to analysis of any component of A(r), at a fixed altitude. That component is 
expressed in terms of spherical harmonics as follows 

90 n • 

Aj(6,<l>) = X X £ c nm cos ( m( l ) ) + S n m sin(m<|>) ] P™ (cos 6) 

n= 1 m = 0 


( 8 ) 


90 


= X 

n=l 


n i i m 

£ [ A nm cos(m<t» + a nm ) ] P„ (cos 9) 

m=0 


where Ai is the i lh component of A(r) at colatitude 0 and east longitude <{>, 

p m (cos0) are the Schmidt normalized associated Legendre functions of 

n v J 

degree n and order m, and C n m and S n m are the spherical harmonic 
coefficients of the magnetic anomalies. A n m = [Cnm 2 + Snm 2 ]^ 2 is an 
amplitude factor and a n m = tan-^-Snm/Cnm)* a phase factor. The 
corresponding power spectrum is defined by 

n p *7 i c 2 

D _ y 1 - nm — ° nm / q \ 

Fn -2* 2n+l k ; 

m=0 

The degree correlation between two such spherical harmonic analyses is 
defined by 


n 


X (Cnm C'nm + S n m S'nm) 
m=0 


n n 

{ [ X (Cnm 2 + S nm 2 )] [ X (C'nm 2 + S’nm 2 )] } ^ 2 
m=0 m=0 


( 10 ) 


Given the analysis of dawn and dusk data according to equation (8), a 
correlation analysis can then be performed to isolate the "common" 
features of the dawn and dusk data, where "common" has an objective 
definition. To define "common", it is noted that, as in pass-by-pass 
correlation, a correlation coefficient for an individual degree and order 
between the two analyses can be defined as 


P am — COS(a n m ■ Ot'nm)* 


( 11 ) 


The criteria chosen to identify "common" features between dawn and 
dusk is that p n m ^ 0.7 and A nm and A' n m not differ by more than a factor 
of two. When this criteria is not met, the terms from that degree and order 


are deleted from the analyses. This criteria was chosen after trying a suite 
of combinations. The choice is admitedly subjective. It is judged to 
preserve most of what are regarded as "common" features while 
succesfully rejecting most of the noise. 


Advanced Processing 

Advanced processing consists of various steps as dictated by physical 
processes controlling D(r,t) in a particular region and our ability to remove 
their magnetic effects. Here lies the heart of the present effort. Equatorial 
and polar regions pose specialized problems in modeling and reduction of 
satellite magnetic data and, therefore, these regions are treated separately. 
A generalized processing sequence (Figure 3) is as follows: an application 
and removal of physical/empirical correction for the appropriate 
ionospheric system (for equatorial electrojet, vector: Langel et al., 1993, 
scalar: Ravat and Hinze, 1993; for polar regions, this study), identification 
and retaining only the correlative components between nearby passes, i.e., 
pass-to-pass correlation, removal of cross-line differences between dawn 
and dusk intersecting profiles (Taylor and Frawley, 1987; Ravat, 1989), 
equivalent source inversion for the purpose of altitude-normalization 
(Dampney, 1969; Mayhew, 1979, von Frese et al., 1981; Langel et al., 1984, 
von Frese et al., 1988), spherical harmonic covariant correlation between 
dawn and dusk Magsat maps. A final equivalent source inversion is then 
used to compute common-altitude anomaly maps and effective volume 
susceptibility constrasts of the lithosphere. 

Separation of fields from ionospheric and lithospheric sources. As noted, 
until recently the removal of D(r,t) has largely been treated as a signal 
processing problem. A(r) was considered to be the consistent part of the 
residual of the data processing sequence. Changing to a technique of 
modeling and removing D(r,t) was pioneered by Yanagisawa and Kono 
(1984, 1985) followed by Cohen and Achache (1990), whose methods have 
been modified and extended for part of the present analysis. The change to 
modeling and removing D(r,t) from data, even in part, requires effort to 


ensure that the ionospheric model is not fitting a significant portion of A(r) 
as well. This issue is important, and difficult, because there is no clear or 
unambiguous separation between ionospheric and lithospheric magnetic 
anomalies in terms of their amplitude or wavelength characteristics. 

In the present study, separation of D(r,t) from AB(r,t) to isolate A(r) 
is carried out in a bootstrap, iterative manner. First an initial estimate of 
the crustal field was removed from AB(r,t). Then, for each pass of retained 
data, an estimate of the form of D(r,t) was made along the track of the 
pass. The amplitude of that form was then determined by least squares fit 
of the form to the pass and the resulting, scaled D(r,t) subtracted from 
AB(r,t) to give an estimate of A(r). 

In the case of equatorial/mid-latitude data, the process utilized 
heavily the modeling of fields from the equatorial electrojet by Ravat and 
Hinze (1993) for scalar data and by Langel et al. (1993) for vector data. 

For vector data, the underlying assumption of Langel et al. (1993) is 
that the sources and fields in question are organized with respect to the 
dip equator (see, e.g.. Figures 1 and 2). First the data are filtered, pass-by- 
pass using a 12000 km cutoff high pass Kaiser filter. Langel et al. (1993) 
examined the extent of equatorial electrojet fields in dusk and dawn data 
and concluded that such effects were absent in dawn data except as dawn 
residuals were contaminated by the field model used to form the residuals. 
Accordingly, the initial estimate of A(r) was derived from Magsat dawn 
data. 

The remaining procedures include fitting data from longitude swaths 
of 45° under the assumption that any remaining effects of A(r), or other 
contamination, will tend to be averaged out. This procedure was followed 
for both dawn and dusk data sets which were then transformed into dip 
latitude coordinates, where dip latitude is defined as 5 = tan-^O.Stanl), 
where I is the geomagnetic inclination evaluated using the GSFC( 12/83) 
field model (Langel and Estes, 1985a). One of two fitting procedures is then 
employed. In the first, the ABe and AB r components are represented by the 
appropriate derivatives of a potential function of the form 


V 


( 12 ) 


1 5 


= a 


N t 

E 

n= 1 


( ”) G n Pn (cos 0d) 


where 0d = 90° - 5. Care was taken to avoid edge-effects at the poleward 
extremes of the data and Nt was chosen after a principal components 
analysis indicated which of the solution coefficients, the G n , were 
significant. The G n were dependent upon the daily sunspot number, Rd, 
according to G n = G' n + G" n *Rd- 

For AB<j> the data were fitted with the function 

1 Nt 

AB(j) = y , I lO'i® X Gn (cos 0d) (13) 

v ; n=0 


where h is the altitude of the data point above the mean radius of Earth, 
and a is an empirically determined constant. Similar attention was given 
to edge-effects and Nt as for the analysis of equation (12). 

A moving window approach in longitude was adopted with 128 

windows, each of 45° width. This results in a shift of 2.8125° between 
windows. Seasonal variation was investigated by dividing the data into 
overlapping subsets of two month duration. Figure 5 shows the ABq and 
AB r calculated from Equation (12) as based on data from March-April, 

1980, in the longitude band 240°-300°. The error bars are the rms of the 

residual data to the field computed from Equation (12). 

To "correct" AB(r,t) for the effects of D(r,t) remaining after filtering, 
the appropriate two month period, sunspot number, and longitude bin 
from the data were identified and the D(r,t) from the analysis of Langel et 
al. (1993) computed along each pass. That D(r,t) was then scaled by 

keeping its form but scaling its amplitude by least squares fit to the 
measured AB(r,t). The resulting scaled D(r,t) was subtracted from AB(r,t). 

A slightly different procedure (following Ravat and Hinze, 1993) was 
used in making ionospheric corrections for scalar equatorial anomaly data. 
The two differences in the procedures are (1) no initial estimate of crustal 
anomalies was removed from the scalar data when modeling ionospheric 
fields, and (2) the longitude swath in averaging scalar data was 90°. As for 


vector data, guide functions of ionospheric field variations are constructed 
in the dip latitude - longitude space using moving window averages. As 
mentioned, the swath width was 90°; the step between windows was 1° 
longitude. As for vector data, the procedure was carried out for dawn and 
dusk data separately. These ionospheric guide functions, which are an 
approximation of D(r,t), were then least-squares scaled to data (i.e., AB or 
AB(r,t)). To avoid over- scaling of D(r,t), in fitting A(r) + D(r,t) in the areas 
of strong lithospheric anomalies, Ravat and Hinze (1993) used an added 
constraint to the linear least-squares problem that the fit to ionospheric 
component must globally reduce cross-over differences between dawn and 
dusk total intensity observations. This method by-passed the effort 
required in first removing an estimate of crustal anomalies from data. 

The above described process was carried out for all components of 
AB(r,t) and for AB(r,t), the scalar residual. 

At high latitudes, i.e., the auroral belt in particular, D(r,t) is much 
more complex than at low and mid latitudes. In particular, AX and AY show 
extremely large and variable behavior in the auroral belts due to the 
presence of field aligned currents. Although an attempt was made, the 
techniques available to us did not prove adequate to isolate what could be 
regarded as a reliable estimate of the anomaly field in these components. 
Accordingly, the high latitude analyses to be described apply only to AB 
and to AZ. Figure 6 illustrates problems faced in the high latitude regions, 
where AB r (i.e., -AZ) component shows distinct anomaly patterns 
systematic in geomagnetic latitudes. It is also apparent from these maps 
that the uncorrected anomalies are also significantly different for the two 
data subsets for each of the polar region. 

Also, in the case of high latitude data, because of its high dependence 
on solar activity and consequent variability, a robust first estimate of 
lithospheric magnetic anomalies could only be obtained through 
quantitative comparisons of data sets having different observational 
characteristics (e.g., different local-times, altitude ranges, etc.). Hence, a 
spherical harmonic covariance correlation of signal from the dawn and 
dusk quiet subsets was used as the initial estimate of crustal anomalies. 

A major difference between the equatorial/mid latitude region 
analysis and the analysis for the polar regions lies in the systematics of the 
auroral ionospheric fields (i.e., their representation in terms of an optimum 


coordinate system). In the case of the equatorial/mid-latitude twilight 
Magsat data, it was possible to organize ionospheric variations with respect 
to magnetic - (dip) latitudes. However, in the polar ionosphere, geomagnetic 
(dipole) latitude, or a refinement thereof, controls the organization of 
electrical currents and consequently their magnetic effects. Another 
important difference, for Magsat data, is that as the satellite orbit 
traversed the polar regions, because it was not in a strictly polar orbit it 
sampled a range of local-times rather than just dawn and dusk. In 
geomagnetic coordinates this difference is even more pronounced because 
of the 11.2° tilt of the geomagnetic coordinate system with respect to 
Earth's rotation axis. Data organization for study and representation of the 
resulting magnetic fields must take into account that ionospheric current 
morphology in auroral and polar regions strongly depends on the magnetic 
local time (MLT). Thus, after subtracting the initial estimate of A(r) from 
AB(r,t), the data averages were compiled in roughly 2° by 2° equal-area 
bins in the dipole latitude - MLT coordinate system. The results for AB for 
the northern hemisphere dawn and dusk data are shown in Figure 7. The 
premise is that such averages are an esti ma te of the average ionospheric 
field, i.e. constitute an D'(r), and can be used to at least partially correct 
AB(r,t) for the effects of D(r,t). 

The correction of AB(r,t) by D'(r) is accomplished similarly to the 
procedure for equatorial and low latitudes. Along each pass the values of 
D'(r) are collected and a common scaling factor, say s, determined for that 
pass such that the scaled D'(r) best fits AB(r,t) in a least squares sense. 
Then the estimate of the crustal field is AH*") = AB(r,t) - D'(r). As 
previously noted, this procedure was carried out for AB and for AZ. Figures 
8 and 9 show ionospheric field corrected dawn and dusk AB anomaly maps, 
respectively (to be compared with Figures 1 and 2). 

Processing to improve signal-to-noise ratio, altitude normalization, and 
gridding. Continuing the description of the flowchart in Figure 3, the stage 
has been reached where as much has been done as is possible to remove 
the effects of D(r,t) from AB(r,t). However, at this stage, there still are 
apparent inconsistencies between nearby passes and between dawn and 
dusk data sets. The premise for further data processing is that lithospheric 
or crustal anomalies should be consistent (not necessarily the same) 



between nearby passes and should be identical (in practice, similar) 
regardless of measurement time. The purpose of this phase of the 
processing is to extract such common components from data (signal) and 
suppress variable components (noise). It is important to note that if the 
ionospheric components have similar wavenumber and amplitude 
characteristics between two data sets being compared, the processing will 
not be able to distinguish their nature: at this stage, it is assumed that this 
distinction was made by the earlier corrections. 

We have compiled maps and statistical parameters showing 
improvements resulting from our processing at various processing steps 

and each of the components (e.g, dawn-dusk difference maps, maps of 
correlation between dawn and dusk 220 km x 220 km equal-area 
averages, and histograms of the above average anomaly differences). 
However, it is impractical to show all of these parameters at each of the 
processing steps and for each component. And hence, only the 

improvements due to ionospheric field corrections are shown in the map 
form for AB 220 km x 220 km equal-area averages; the parameter chosen 
here is the dawn-dusk difference (Figure 10a and b). The remaining 

processing improvements are shown, once again for AB, as histograms of 
the above dawn-dusk average differences (Figure 11). 

To improve the consistency of data (signal) between two nearby 
passes, the pass-by-pass correlation technique (Alsdorf, 1991), described 

earlier, was used. As noted, the commonalities between two equal- length 

profiles can be quantitatively compared, i.e., a correlation coefficient 

derived, for each wavenumber. Admittance criteria for a particular degree 
of correlation, however, are qualitative. For example, two profiles with 
identical spatial locations should be required to correlate with high degree 
of correlation; as the distance between the profiles increases, the 

assumption of two dimensionality of anomalies perpendicular to the 
direction of the profiles is lost and the wavenumber with lower correlation 
coefficients cannot be automatically rejected as noise. After 
experimentation with a range of correlation coefficient cut-off criteria and 
their effect on all the different segments of the anomaly maps, 

wavenumbers with correlation coefficients of 0.3 and above were retained 
for all dawn/dusk and equatorial/polar data segments. A considerable 
improvement has been realised due to the implementation of the pass-by- 


pass correlation technique for each of the data subsets (Figures 12 and 13; 
to be compared.with Figures 8 and 9) 

After completing adjustments to the data from intercomparing 
adjacent passes, an analysis was conducted of the crossover differences 
between dawn and dusk passes. The purpose of adjusting cross-line 
differences (Taylor and Frawley, 1987; Ravat, 1989) was to adjust 
simultaneously the levelling inaccuracies between the profiles within the 
individual data sets and between the dawn and the dusk subsets. In this 
study, the cross-line adjustments were implemented only on 
equatorial/mid-latitude data. 

At this point the, multiply corrected, data are spaced unevenly in 
position, including altitude. This introduces artificial noise when the data 
are averaged and when comparing dawn and dusk data which are at 
differing altitudes. Also, processing of the three components at each local 
time has to this point proceeded independently. The next step in the 
process is to derive equivalent point dipole representations which can be 
used to reduce the data to common elevation and to facilitate further 
intercomparison of data from the two local times. 

Derivation of this representation was inspired by Dampney's (1969) 
method for synthesizing Bouguer gravity measurements on an irregular 
three dimensional grid. The synthesis consisted of a mathematical 
representation of the data in terms of discrete point masses at some, 
arbitrary, fixed depth below the Earth's surface. Mayhew (1979) adapted 
this method to the synthesis of magnetic anomaly data acquired by the 
POGO satellites. In this method the satellite magnetic anomaly data are 
represented by an array of dipoles at the Earth’s surface. The dipoles are 

assumed to be aligned along the direction of the Earth's main field, as 

determined by a spherical harmonic model, and their magnitudes are 

determined so as to best reproduce the anomaly data in a least-squares 
sense. Following Dampney, this is called an equivalent source model. 
Working independently, von Frese et al. (1981) proposed the same method 
of analysis, only applied to gravity as well as magnetic anomaly data. 

Such inversions were performed, separately for dawn and dusk data, 
for the combined, corrected, AZ and AX data, except at high latitudes where 
only AZ were included. Inversion was also performed for the AB data. The 
AY component was not used for two reasons. First, because even after 
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application of the corrections described above, AY for dusk was judged to 
still be contaminated to an unacceptable amount. Second, because the very 
clean dawn AY, at low and mid latitudes, can be used as an independent 
test of the validity of the equivalent source solution. Examples of AB dawn 
and dusk anomalies at 400 km elevation after the inversion are shown in 
Figures 14 and 15, respectively. 

The inversion also allows recomputation of anomalies at common 
altitude and at specified grid intervals as required for the application of 
Spherical harmonic covariant correlation previously described. Spherical 
harmonic covariant correlation provides a powerful tool for further 
comparison and consistency improvement between the dawn and the dusk 
data sets. As mentioned earlier, the technique is used for retaining 

commonalities in maps rather than profiles and suppressing the 

contribution made by undesirable harmonics (as apparent from their small 
correlation coefficients and differences in amplitudes). Figure 16 shows 
correlation coefficients between the dawn and the dusk anomaly maps for 

each harmonic degree for each of the anomaly components before and 

after the selection of criteria for common characteristics (correlation 
coefficients £ 0.7 and amplitude factors A nm 's within factor of two). As 
described earlier, coefficients that do not meet these criteria are deleted 
from the analysis. It is apparent from the spectra in Figure 16 that before 
the above selection the correlation coefficients are much lower than 1.0 
(solid lines), indicating that there is residual external field contamination 
in the two maps. After the selection, however, the two signals (maps) are 
seen to correlate to a high level (long-short dashed line in Figure 16). 
Figures 17 and 18 show the power at each degree harmonic for dawn and 
dusk maps, respectively, and for each of the components before and after 
the process. It can be observed from these figures that most of the original 
power is retained by the analysis, reduction of power being limited to non- 
correlative parts between the two corresponding data sets. 

After the spherical harmonic correlation analysis, all processed data 
sets were subject to simultaneous inversions (vector and scalar separately) 
from which final vector and scalar anomalies have been computed at 400 
km (Figure 19). 


SUMMARY 


Vector and scalar magnetic anomaly maps have been prepared from 
Magsat data using a systematic approach that integrates a variety of useful 
processing schemes developed during the last decade. The approach 
initially focuses on modeling and removal of ionospheric magnetic fields in 
both equatorial and polar regions and subsequently relies on signal 
processing techniques to minimize remaining discrepancies between data 
measured at different local-times, the premise being the lithospheric 
anomalies should be consistent regardless of the measurement time. This 
technique has resulted in identifying the most robust Magsat magnetic 
anomalies that are most likely lithospheric in origin (Figure 19). As a 
corollary, the anomalies that do not appear on these maps but appear in 
other published maps that have been prepared without taking into account 
the ionospheric contamination are less robust and therefore such 
anomalies should be carefully scrutinized for ionospheric contamination 
and repeatability before they are interpreted. 

In addition to the mapping of the most robust satellite magnetic 
anomalies, new techniques have been devised to model polar ionospheric 
fields in AZ and AB components and to map polar magnetic anomalies from 
satellites that were hitherto untenable because of the large ionospheric 
contamination in the polar regions. Moreover, for the first time, analytical 
vector equatorial ionospheric corrections were used to map equatorial 
lithospheric magnetic anomalies. 

Thus, significant progress has been made during the award period in 
analyzing and modeling ionospheric fields in Magsat data and devising 
techniques to map most repetitive Magsat magnetic anomalies of the 
lithosphere. The techniques developed during this project are also 
applicable for the analysis of future low altitude satellite magnetic data 
sets and would prove invaluable in their initial analysis. 
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Table 1. Processing parameters (see processing flowchart. Fig. 3): 


Process 

Equatorial vector 

Equatorial scalar 

Polar 

Quiet selection 

Kp £ 1 + 

Kp < 2+ 

AE < 50 nT 
ABvar < 80 nT 2 

Main Field 

Dawn 6/83-6 

Dawn 6/83-6 

GSFC 12/83 

High-pass (^.) 

12000 km 

4000 km 

4000 km 

No. of Iterations 




of ionospheric 
field removal 

1 

N/A 

1 

Pass-to-pass 

cc £ 0.3 

N/A 

cc > 0.3 

Cross-line 

adjustments 

zero order 

zero order 

none 

Covariant spher. 
harmonics 

cc > 0.3 

Ampl. fact. = 2 

cc ^ 0.3 

Ampl. fact. = 2 

cc £ 0.3 

Ampl. fact. = 2 
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Figure 1 

Figure 2 

Figure 3 
Figure 4 

Figure 5. 

Figure 6. 
Figure 7. 

Figure 8. 
Figure 9. 


FIGURES 


Total intensity anomaly (residual) map from Magsat dawn data 
prior to ionospheric corrections. Units are nT. Cl = 2 nT. Van 
Der Grinten projection. 

Total intensity anomaly (residual) map from Magsat dusk data 
prior to ionospheric corrections. Units are nT. Cl = 2 nT. Van 
Der Grinten projection. 

Flow chart of Magsat data processing. 

Illustration of correction of attitude jump in Magsat data, x x x 
indicates original AZ data; solid line indicates corrected data. 

Computed values of AZ (top) and AX (bottom) from Equation 
(12) for the longitude range 240-300° during March- April, 
1980. The ordinate is dip latitude; the abcissa is in nT. The 
error bars reflect the 1 standard deviation scatter of the data 
used in the fit from the values computed from Equation (12). 

AB r component anomalies over polar regions showing the 
bands of anomalies systematically along geomagnetic latitudes, 
(a) north polar, dawn; (b) north polar, dusk; (c) south polar, 
dawn; (d) south polar, dusk. 

Averages of ionospheric field, D'(r,t) = AB(r,t) - A'(r), for AB 
component in MLT-geomagnetic (dipole) latitude coordinates. 
A’(r) is the initial estimate of A(r). Units are nT. Cl = 2 nT. (a) 
north polar, dawn; (b) north polar, dusk; (c) south polar, dawn; 
(d) south polar, dusk. 

Ionospheric field corrected anomaly map of AB from Magsat 
dawn data. Units are nT. Cl = 2 nT. Van Der Grinten projection. 

Ionospheric field corrected anomaly map of AB from Magsat 
dusk data. Units are nT. Cl = 2 nT. Van Der Grinten projection. 
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Figure 10. Difference map of dawn-dusk 220 km x 220 km equal-area 
averages of AB. Cl = 2 nT. Van Der Grinten projection, (a) before 
ionospheric corrections; (b) after ionospheric corrections. 

Figure 11. Histograms of AB dawn-dusk differences before and after 
ionospheric corrections, (a) north polar; (b) equatorial; (c) south 
polar. 

Figure 12. Dawn AB anomaly map after pass-by-pass correlation. Cl = 2 
nT. Van Der Grinten projection. 

Figure 13. Dusk AB anomaly map after pass-by-pass correlation. Cl = 2 nT. 
Van Der Grinten projection. 

Figure 14. Dawn AB anomaly map at 400 km after equivalent source 
inversion. Cl = 2 nT. Van Der Grinten projection. 

Figure 15. Dusk AB anomaly map at 400 km after equivalent source 
inversion. Cl = 2 nT. Van Der Grinten projection. 

Figure 16. Degree correlation between final dawn and dusk maps, (a) AZ; 
(b) AX; (c) AY; (d) AB. 

Figure 17. Power spectra of dawn anomaly maps, (a) AZ; (b) AX; (c) AY; (d) 
AB. 

Figure 18. Power spectra of dusk anomaly maps, (a) AZ; (b) AX; (c) AY; (d) 
AB. 

Figure 19. Combined corrected anomaly map from Magsat data. Units are 
nT. Cl = 2 nT. (a) AZ; (b) AX; (c) AY; (d) AB. 
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FLOW-CHART OF MAGSAT DATA PROCESSING 

(■ : keyed lo parameters In Table I) 
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